function [ct,pt] = blstheta(so,x,r,t,sig,q) 
%BLSTHETA Black-Scholes sensitivity to time until maturity change. 
%   [CT,PT] = BLSTHETA(SO,X,R,T,SIG,Q) returns sensitivity in option value with
%   respect to time.  SO is the current stock price, X is the exercise price,
%   R is the risk-free interest rate, T is the time to maturity of the option
%   in years, SIG is the standard deviation of the annualized continuously
%   compounded rate of return of the stock (also known as volatility), and
%   Q is the dividend rate.  The default Q is 0.  CT is the theta of a call
%   option, and PT is the theta of a put option. 
%       
%   Note: 
%     This function uses normpdf, the normal probability density function 
%     and normcdf, the normal cumulative distribution function in the
%     Statistics Toolbox. 
% 
%   For example, [c,p] = blstheta(50,50,.12,.25,.3,0) returns  
%   c = -8.9630 and p = -3.1404. 
% 
%   See also BLSPRICE, BLSDELTA, BLSGAMMA, BLSRHO, BLSVEGA, BLSLAMBDA. 
% 
%   Reference: Hull, Options, Futures, and Other Derivative Securities,  
%              2nd Edition, Chapter 13. 
 
%       Copyright 1995-2006 The MathWorks, Inc.
%       $Revision: 1.6.2.4 $   $Date: 2009/04/15 23:07:23 $ 
 
 
if nargin < 5 
  error('finance:blstheta:missingInputs','Missing one of SO, X, R, T, and SIG.') 
end 
if any(so <= 0 | x <= 0 | r < 0 | t <=0 | sig < 0) 
  error('finance:blstheta:invalidInputs',sprintf('Enter SO, X, and T > 0. Enter R and S >= 0.')) 
end 
if nargin < 6 
 q = zeros(size(so)); 
end 
 
message = blscheck('blstheta', so, x, r, t, sig, q);
error(message);


% Perform scalar expansion & guarantee conforming arrays.
try
    [so, x, r, t, sig, q] = finargsz('scalar', so, x, r, t, sig, q);
catch
    error('Finance:blstheta:InconsistentDimensions', ...
        'Inputs must be scalars or conforming matrices.')
end

% blspriceeng works with columns. Get sizes, turn to columns, run engine,
% and finally turn to arrays again:
[m, n] = size(so);

% Double up on fcn calls since blsprice calculates both calls and puts. Do
% this only if nargout>1
NumOpt = numel(so);
callSpec = {'call'};
callSpec = callSpec(ones(NumOpt,1));
putSpec = {};
if(nargout)>1
    putSpec = {'put'};
    putSpec = putSpec(ones(NumOpt,1));
    % double up the rest of the input args
    [so, x, r, t, sig, q] = deal([so(:);so(:)], [x(:);x(:)], [r(:);r(:)], ...
         [t(:);t(:)], [sig(:);sig(:)], [q(:);q(:)]);
end

OptSpec = [callSpec;putSpec];
OutSpec = {'theta'};

% call eng fuction
theta = blspriceeng(OutSpec, OptSpec, so, x, r, t, sig, q);

% Now separate calls from puts
ct=reshape(theta{1}(1:NumOpt), m, n);
if(nargout>1)
    pt = reshape(theta{1}(NumOpt+1:end), m, n);
end